Free energy fluctuations in Ising spin glasses 
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The sample-to-sample fluctuations of the free energy in finite-dimensional Ising spin glasses are 
calculated, using the replica method, from higher order terms in the replica number n. It is shown 
that the Parisi symmetry breaking scheme does not give the correct answers for these higher order 
terms. A modified symmetry breaking scheme with the same stability is shown to resolve the 
problem. 



PACS numbers: 75.50.Lk, 05.50.+q 

The Parisi replica symmetry breaking scheme jj], |^, 
H was a milestone in the study of spin glasses. It was 
proposed as long ago as 1979 and has been extensively 
used ever since. It is the purpose of this Letter to point 
out that an important modification of it is needed in 
order to calculate the sample-to-sample fluctuations of 
spin glasses. 

The replica method is usually used in spin glasses to 
calculate the average free energy F as a function of in- 
verse temperature p from the partition function Z via 



f3F 



-InZ 
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The overbar means averaging over bond configurations. 
Apart from the free energy, the replica method also gives 
in principle access to other physical quantities as well. 
Expanding the logarithm on the right hand side in Eq. (^) 
one gets 



lnZ n = —n/3F + — i'AF 2 + 



(2) 



where AF = (In Z — InZ )//? is the mean-square 
sample-to-sample fluctuation of the free energy, and the 
coefficients of the higher order terms are higher order 
cumulants. In order to obtain the coefficients in this ex- 
pansion from the replica method, it is necessary to keep 
the replica number n small but finite throughout the cal- 
culation (as opposed to the case of the free energy itself 
where it is possible to set n = early on), which makes 
it rather cumbersome. Moreover, using Parisi's replica 
symmetry breaking scheme ||, ||, |[ Q], it will be shown 
below that the coefficients so obtained differ when as- 
suming n positive or negative. Note that the replica 
method intrinsically requires n to be non-integer, and 
that it is no more unnatural to take n negative than it 
is to take it non-integer. Thus conflicting answers are 
obtained for the same physical quantity which indicates 
that there is a problem with Parisi's replica symmetry 
breaking scheme. We will argue that the correct answer 
is in fact given by the n < solution. Furthermore we 
show how the symmetry breaking scheme can be modified 
to give the correct results when n > 0. However, this new 



replica symmetry breaking scheme gives the same results 
as Parisi's scheme for quantities like the free energy or 
the distribution of spin overlaps P(q), and has the same 
stability. 

We start from the usual replica field theory for d- 
dimensional spin glasses as derived, e.g., in j^] or [[|, 
where we have a free energy functional 
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and 



expi-HrepUap})- (4) 



Here, q aa — 0, we have omitted some unimportant terms 
of order g 4 , and set r = 1 — T/T c . The fourth order 
term included is the one responsible for replica symmetry 
breaking. We assume that the dimension d is above the 
special dimension 8 to keep the loop expansion straight- 
forward l|. 

Since we are keeping n finite throughout, it will be 
necessary to review some of the results that have been 
derived in § @, |, | for into ger n, allowing us to 
write, e.g., q a p with integer indices a and /3, and of results 
that have been derived for arbitrary (noninteger) n and 
with an infinite number of symmetry breaking steps, such 
that we must re- write q a p in terms of q(x) where x is 
a continuous variable between n and 1. We will have 
to switch between these notations frequently but note 
that they are equivalent. In order to distinguish better 
between the cases n > and n < 0, we will from now on 
use the notation n + for positive n and n~ for negative 
n. The letter n itself will be used when no distinction is 
needed. 
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FIG. 1: The two types of g-function for positive (n + ) and 
negative (n~) replica number. For x > 3n + /2, the two curves 
coincide. 



To Gaussian order we get from Eq. 



In Z n = 



-7~Lrep{qa^} 
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(5) 



where is the Parisi saddle point solution, k is a d- 
dimensional wave vector and A^d^ are the eigenvalues 
of the Hessian, evaluated at the saddle point solution, 
and their degeneracies. 

The first (mean field) term in Eq. (||) has been worked 
out by Kondor Q for finite n + . We repeat his calculation 
here for n~ since it gives a qualitatively different result 
which we will need later on. The free energy at mean- 
field level for the Parisi function q n {x) at finite n is (cf. 
Eq. (9) Of §) 



~ dxql(x) + j^j dxq*(x)- 



-J dx \xq 3 n (x)+3q 2 n (x)J dtq n (t)\ , (6) 

when extremized with respect to q n (x) . The solutions for 
n + and n~ are 

3wn+/4y n+<x<3n+/2 
q n +(x) — ^ wx/2y 3n + /2 < x < X\ (7) 

wx\/2y x\ < x < 1 

n~ < x < 

q n -(x) = { wx/2y < x < xi , (8) 
wx\/2y x\ < x < 1 



where x\ = 1 — \/l — Ayr/w 2 is the usual breakpoint of 
the Parisi g-function. Fig. [I] shows the two solutions for 
illustration. 

We note that the kind of problem we are going to en- 
counter already shows up at this level: Ti rcp [q n ~ (x)] has 
no terms of higher order than n~ , while Ti. lcp [q n + (x)] 
does have additional terms of order n +6 and higher M. 



We also note that at mean-field level, the variance AF 2 
is since there is no n 2 -term. 



In order to calculate the fluctuation corrections, 
useful to define 



/ = ^^ln(fc 2 + A p 



it is 



(9) 



By differentiating / with respect to k 2 one obtains 
dl ^ d 



d{k 2 ) 



E 



fc 2 + A t 



(10) 

Q</3 



where G a p !7 s are the propagators, which are essentially 
the inverse of the Hessian. The propagators have been 
calculated exactly by De Dominicis et al. || in the "con- 
tinuum limit" of infinitely many replica symmetry break- 
ing steps. While their results are too long to be quoted 
here, we mention that the diagonal propagators G Q /3, a /3 
are in this limit denoted by Gff , are labelled by a con- 
tinuous variable x € [n, 1], and we have 



- j dxGfi. 



(11) 



For ease of notation, we will from now on drop the sub- 
script 11 from the propagators (which is superfluous for 
our purposes) and replace it by n to indicate that the 
propagators here still depend on it. While the propaga- 
tors in H were derived for n = 0, it is possible but tedious 
to extend the calculation to finite n. Fortunately, there is 
a simple argument to obtain the propagators's exact form 
for n~ . Since they are labelled by the variable x, which 
is nothing but the inverse of the Parisi function, i.e. x(q), 
they are effectively labelled by q. From this we can infer 
that G^L = const, for those x where q n - (x) is a con- 
stant, therefore (cf. Eq. (§)) = G°° for x < 0, and 
for x > 0, G^L = Gq x because q n - (x) = go(^)- We can 
thus express the propagators for n~ entirely in terms of 
the n = propagators. For n + this simple argument does 
not work because q n + (x) from Eq. (Q) is (almost) iden- 
tical to the g-function in a field, and it has been shown 
in H that this leads to an additional n + -dependent term 
in the propagators. We will see below that luckily we do 
not need to work this out in detail. 

We are now in a position to calculate AF 2 for rT . 
First, we calculate 



dl 
<9(fc 2 ) 



dxG x * 



dx Gq x 



on 



-G, 
(12) 



Thus there are no terms of higher order than n~ . This 
already implies, by comparison with Eq. (^), that the free 
energy fluctuations have a Gaussian distribution since all 
higher order cumulants are zero. In order to calculate the 
variance of this distribution, A_F 2 , we only need to know 
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the coefficient of the n -term, which is simply Gq°, and 
integrate it with respect to k 2 . 
given by §] 



The propagator Gq is 



1 dt d 2 
Tds~di 
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f(s,xi) + — — f{xi,t) ) + f(xi,xi), 



where 



f(s,t) 



k 2 + yq 2 (s) + yq 2 {t) 



(13) 



(14) 



is the inverse of the x — eigenvalues of the Hessian 
from the replicon sector The function / can easily 
be integrated with respect to k 2 jl3| , and the integrals in 
Gq° can then be worked out, resulting in 



mi 



J := / d{k 2 )G^ 
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(16) 



^yk 2 \J Aqjk 2 + w 2 x\ \f^yk 2 + w 2 x\ ' 



such that 



^2 AF 2 = _Y_ 



d d k 
(2n) d 



(17) 



The integrals over k diverge unless wc introduce a cutoff, 
which is tacitly implied in Eq. (|l7|). 

According to what we found so far for n~ , the free 
energy fluctuates with a Gaussian distribution and a 
variance as given by Eq. @. This is, we believe, the 
physically sensible solution: if, by the usual argument, 
a spin glass sample is divided into many subsystems, 
each should contribute a random dominant bulk term 
and a subdominant surface term to the free energy, and 
by the central limit theorem this should give rise to a 
self-averaging, Gaussian-distributed quantity. 

If we repeat the above procedure for n + , however, we 
find the following situation, which is illustrated most eas- 
ily in a large-fe 2 expansion of the propagators || , 



G 



a/3, a/3 



1 2r + 2yql 
k 2 fc 4 



+ 



((2r + 2yq 2 afj ) 2 + w 2 ((q 2 ) aa + (q 2 )p - 2q 2 a0 )) 

+ 0(l/k s ). (18) 



This expansion allows for a simple evaluation of 
Ea<^^.^) or rather -n/2 f^dxG* x , term by term. 



While the first two terms are identical for n + and n 
the coefficient of the l/fc 6 -term is 



-(4T 2 -8yrg+ ^(1-4x1/5)) 



^(Ar 2 -2w 2 q)-n+ 6 



8lw 4 



5 ■ 128y 2 ' 



(19) 



while for n one obtains the same but without a 
ri~ -term. Here we have used the abbreviation q = 
~ In 1n( x ) dx, which is independent of n. This shows 
that at least one expansion coefficient, and thus the phys- 
ical consequences, changes when changing the sign of n, 
since a nonzero n 6 -term implies a non-gaussian probabil- 
ity distribution of the free energy. 

The difference between the two cases n + and n~ can 
be eliminated by modifying Parisi's original symmetry 
breaking scheme for n + in the following way. The q- 
matrix is divided into p boxes on the diagonal, each of 
which contains a Parisi-type symmetry broken matrix 
Qa/3, while the rest of the matrix is zero. This is il- 
lustrated in Eq. (p0|) for p = 4, 



Qa/3 
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Qa/3 



\ 



Qa/3 



(20) 



Now we first let the number of symmetry breaking steps 
R go to infinity, then we let p go to infinity, all the while 
keeping n + finite. 

This procedure is in fact very closely related to, 
but certainly not identical to, the original Parisi 
scheme. If we introduce the usual "integers" n + = 
mo, mi, . . . , rtiR, rtiR+i = 1 which characterise the block 
sizes of the symmetry broken matrix, then the differ- 
ence between this scheme and the original one is that 
mi = n + /p is allowed to go to before mo- This way 
of looking at it shows precisely where the difference be- 
tween positive and negative n lies in the original scheme: 
for n~ , mi does go to before mo since now lies in the 
interval 1]. It also shows that we may still use those 
exact results of [pi |^] for the modified scheme which were 
derived for arbitrary m's. 

From Eq. (j^) it follows now that on mean-field level, 
T^rcp [q n + ( x )] is to be replaced by pWrep[Sn+ l v {x)] (there 
are p Parisi blocks of size n + jp each), and in the limit 
p — > oo this kills the n + and higher order terms in 
Kondor's solution (terms of order n + are replaced by 
p(n + /p) m — > (p — > oo) for m > 1), while only the 
term linear in n + survives, unaltered. Thus this proce- 
dure has cured the inconsistency on the mean-field level 
without affecting the average free energy itself. It also 
cures the discrepancy on the Gaussian level, as the fol- 
lowing argument shows. The eigenvalue equation for the 
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Hessian which follows from Eq. (g) is 
(A + 2r) f a p + w qp 7 / 7a +w + 2yq* f a p = 0, 

7 



5 



where A is the eigenvalue and f a a is the eigenvector. 
Since q a /j is zero in the off-diagonal blocks, the eigenvec- 
tors can be chosen to be zero everywhere except in one 
block (and its counterpart on the opposite side of the 
diagonal, if the block does not happen to be on the diag- 
onal). If f a p is nonzero in a block on the diagonal, the 
problem reduces to the eigenvalue problem of the origi- 
nal Parisi matrix, but of size n + /p. There are p blocks 
on the diagonal, i.e. we have to replace — \~ I + dxG^ 
in Eq. (|ll|) by ~p^ f* + /p dx /p , which again in the 

limit p ~ + oo kills all terms of order n +2 and higher 
but leaves the linear term unchanged. There are how- 
ever additional eigenvectors which are nonzero in one of 
the off-diagonal blocks, of which there are p{p — l)/2. 
According to |j, |8| the corresponding eigenvalues are 
again just the x — eigenvalues from the replicon sec- 
tor, given by l//(s,<) from Eq. (|l4|), and their contri- 
bution to dl/d(k 2 ) can once more be reduced to the 
propagator G™ and is (p(p - l)/2)(n+ /p) 2 G° n ° + /p -» 

(n +2 /2)Gq° (p — * oo), while higher order terms vanish. 
Collecting the pieces we get 



dl 

d{k 2 ] 



dxG% x 



-Lrn 



(22) 



Comparison with Eq. (12) shows that this modified type 
of symmetry breaking gives identical results to the ones 
derived before for and thus eliminates the discrepan- 
cies between n + and n~ limits. 

The argument above also shows that the stability of 
the solution is not affected by the modification since all 
eigenvalues of the Hessian are the same as for the original 
scheme and are thus > 0. 

We conclude with a few remarks. A consequence of our 



modification of the symmetry breaking scheme is that 
the Parisi g-function is not well defined for n > since 
the numbers mo, mi, ... do not form a monotonous se- 
quence any more. Instead, it turns into the multi-valued 
function sketched in Fig. ^|, elegantly restoring symmetry 
with q n - (x), which was violated in Fig. |l|. Integrals over 
functions g of q n + (x) have then to be interpreted as 

dxg(q n +(x))= dxg(q l n+ (x)) + / dx g{q\ + {x)), 

Jn+ JO 

(23) 

where the arrows denote the branch in an obvious way. 
Physical quantities like the probability distribution of the 
order parameter, P(q) = dx/dq, are unaffected by this 
manoeuvre since they are only meaningful at n = 0. 

The need to modify Parisi's original symmetry break- 
ing to the scheme described here probably escaped no- 
tice as the two schemes become identical in the limiting 
case n = 0, which is usually considered. Indeed higher 
order terms in the expansion in n have rarely been in- 
vestigated for spin glasses as they are extremely difficult 
to compute when using the established rt + -formalism as 
in |(|. With our modified scheme, however, (or, equiv- 
alently, by working at n < 0) not only are the higher 
order terms guaranteed to be correct but as shown here 
they also become relatively accessible. This may be im- 
portant for future applications, e.g. for computing inter- 
face energies in spin glasses (see O] for an attempt in 
this direction). We would also expect that the modifica- 
tion of the Parisi scheme described here is needed for the 
Sherrington-Kirkpatrick infinite range model. Whether 
this is the case or not could be investigated from the sam- 
ple to sample fluctuations of the ground state energy of 
this model, which can be determined numerically when 
the number of spins N is small |3 . 
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FIG. 2: The multi- valued function q n + (x) for the modified 
symmetry breaking scheme. Note the similarity with q n - (x) 
from Fig. [I]. 
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